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ABSTRACT 


The  object  of  this  work  was  to  explore  the  possibility 
of  scaling  a  Finite  Element  program  to  fit  in  a  small 
desktop  computer.  A  program  originally  written  in  FORTRAN 
IV  and  called  STAP  was  selected  for  the  present  project. 

This  program  is  contained  in  Bathe  and  Wilson's  Numerical 
Methods  in  Finite  Element  Analysis  and  was  conceived  to  solve 
static  linear  elastic  analysis. 

The  program  was  rewritten  in  an  enhanced  form  of  BASIC 
language  suitable  for  the  Hewlett-Packard  System  45  desktop 
computer.  Presently,  only  the  truss  element  is  available; 
however,  the  program  is  structured  so  that  other  element 
types  may  be  easily  implemented.  A  plot  of  the  input  mesh 
is  included  in  the  program  output. 
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I .  INTRODUCTION 


There  is  no  doubt  that  the  finite  element  method 
represents  a  rare  breakthrough  in  engineering  analysis. 

This  technique,  in  the  last  twenty  years,  has  been  sharpened 
and  refined  and  is  today  commonly  utilized  in  stress  analysis 
and  solid  mechanics.  Many  computer  programs  have  been 
written  for  finite  element  analysis  of  static  and  dynamic 
problems,  with  most  of  the  general  purpose  codes  written  in 
FORTRAN  for  large  institutional  computers.  The  last  few 
years  have  seen  the  advent  of  the  personal  computer,  and 
their  memory,  though  limited  at  first,  have  been  rapidly 
increasing  in  capacity.  With  today's  desktop  computers  it 
now  seems  feasible  to  write  a  limited  purpose  computer 
program  which  can  be  implemented  on  a  desktop  computer  that 
has  enough  memory  capacity  to  handle  the  large  memory 
requirements  that  finite  element  analysis  requires. 

This  project's  purpose  is  to  write  a  computer  program 
that,  although  simplified  in  various  areas,  shows  all  the 
important  features  of  more  general. codes.  The 
Static  Structural  Analysis  Program  (SSAP-NPS)  is  a  simple 
computer  program  for  linear  elastic  finite  element  analysis. 
It  is  written  in  an  enhanced  form  of  BASIC  language  and 
has  been  tailored  to  be  utilized  on  a  Hewlett-Packard 
System  45  desktop  computer.  The  main  objective  of  the 
program  is  to  show  that  a  finite  element  similar  in 
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procedures  and  overall  flow  to  more  general  codes 
[References  1,2]  can  be  practical  for  use  on  a  desktop 
computer.  Only  a  truss  element  has  been  made  available  in 
SSAP-NPS;  however,  the  code  can  be  used  for  two-  and 
three-dimensional  truss  analysis  and  additional  elements 
can  be  added  with  general  ease.  The  program  also  has  a 
mesh  plotting  capability  for  verifying  the  input  mesh. 

In  writing  the  SSAP  code,  the  organizational  structure 
and  procedures  of  a  FORTRAN  code  called  STAP  (Static 
Analysis  Program) ,  contained  in  Reference  3,  was  closely 
followed. 
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II.  THE  HEWLETT-PACKARD  SYSTEM  45 
DESKTOP  COMPUTER 

A.  USER  FAMILIARIZATION 

A  photograph  of  the  H-P  System  45  desktop  computer  is 
shown  on  Figure  1.  Users  of  SSAP-NPS  must  have  famil¬ 
iarity  with  its  features,  operation,  and  some  knowledge 
of  its  programming  language,  BASIC  (Beginners  All-Purpose 
Symbolic  Instruction  Code).  References  4,  5  and  6  are 
manuals  which  were  instrumental  in  the  creation  of  the 
SSAP-NPS  program.  Reference  4  contains  all  the  informa¬ 
tion  a  user  would  need  to  utilize  SSAP-NPS.  In  general,  a 
user  must  be  able  to  do  the  following: 

1.  Load  and  run  a  program  stored  on  tape. 

2.  Create  and  size  data  files  on  tape. 

3.  Store  and  purge  files  on  tape. 

4.  Create  and  store  a  program  on  tape. 

B.  FEATURES 

The  H-P  System  45  computer  for  which  SSAP-NPS  was 
created  has  the  following  optional  features: 

1.  Thermal  Printer 

✓ 

2.  Option  203  -  64K  memory 

3.  Optional  tape  drive 

4.  ROM  graphics  package. 

All  of  these  features  are  required  in  using  SSAP-NPS. 


The  requirement  for  64K  memory  may  not  be  absolutely 
necessary  but  since  the  program  itself  requires  about  2 OK 
bytes  of  memory  and  the  only  other  memory  option  larger 
than  20K  is  32K,  it  would  leave  only  12K  of  memory 
available  for  execution  of  the  program.  This  would  limit 
the  program  to  relatively  small  problems. 


III.  GENERAL  FORMULATION  AND  IMPLEMENTATION 


i 
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The  basic  idea  of  the  finite  element  method  is  to 
divide  a  complex  problem  into  a  series  of  simple  inter¬ 
related  problems  that  are  easier  to  solve.  The  whole  is 
modeled  as  an  assemblage  of  discrete  part,  or  finite 
elements. 

In  the  case  of  a  truss  structure,  the  model  is  an 
assembly  of  pin-connected  rods  having  force  and  displacement 
characteristics  known  from  strength  of  materials.  If  each 
interconnecting  member  is  identified  as  an  element,  then 
equations  can  be  written  expressing  the  element  properties. 
The  unknown  forces  and  deflections  can  then  be  found  for 
the  overall  system  by  combining  these  individual  equations 
according  to  the  laws  of  equilibrium,  and  then  solving  the 
resulting  system  of  equations. 

Using  matrix  notation,  the  model  for  a  linear  static 
system  is  expressed  as 

[K]  {U}  =  {R} 

where  [K]  is  the  square  symmetric  stiffness  matrix  of  the 
assembled  structure;  {U}  is  the  column  matrix  of  unknown 
nodal  force  displacements?  and  (R>  is  the  column  matrix  of 
forces  applied  at  the  nodal  points . 
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This  relationship  gives  the  linear  correlation  between 
nodal  forces  and  nodal  displacement  in  the  form  of  a  set  of 
simultaneous  linear  equations  which  is  common  to  problems 
in  static  analysis. 

After  a  structure  has  been  idealized  as  by  above 
relationships,  the  stress  analysis  can  proceed  in  the 
following  three  phases: 

1.  Calculation  of  structure  matrices  [K]  and  {R} 

2.  Solution  of  equilibrium  equations 

3.  Evaluation  of  element  stresses 

The  procedure  for  accomplishing  the  above  analysis  phases 
is  described  only  in  general  terms  in  the  following  discus¬ 
sions.  For  more  complete  discussions  with  illustrative 
examples,  Reference  3  (Chapters  3,  6,  7)  and  Reference  7 
(Chapter  24)  should  be  consulted. 

A.  CALCULATION  OF  STRUCTURE  MATRICES 

The  calculation  of  the  structure  matrices  [K]  and  {R} 
is  performed  as  follows: 

1.  Geometric,  material,  and  loading  data  are  read  and/or 
generated. 

2.  The  element  stiffness  matrices  are  calculated. 

3.  The  structure  matrices  are  assembled. 

1.  Geometric,  Material  and  Loading  Data 

Corresponding  to  each  nodal  point  are  a  number  of 
degrees  of  freedom,  i.e.,  for  the  truss  element  there  are 
three  translational  degrees  of  freedom  at  each  node.  An 
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array  named  the  ID  array  identifies  which  of  these  degrees 
of  freedom  are  actually  used  in  the  analysis.  The  ID  array 
is  of  dimension  number  of  degrees  of  freedom  per  node  *  NUMNP, 
the  total  number  of  nodal  points  in  the  system.  If  ID(I,J)  =0 
the  degree  of  freedom  corresponding  to  the  Ith  degree  of 
freedom  at  nodal  point  J  is  active;  and  if  ID(I,J)  =  1, 
the  degree  of  freedom  is  nonactive.  After  all  active 
degrees  of  freedom  have  been  identified  by  zeros  in  the  ID 
array,  they  are  replaced  by  equation  numbers  by  columns, 
starting  from  1.  The  total  number  of  equations  is  deter¬ 
mined  by  the  maximum  value  in  the  ID  array  and  is  called 
NEQ.  Nodal  location  with  respect  to  a  global  X,Y,Z  set  of 
axes  and  the  loading  at  each  node  are  also  input  at  this 
stage.  In  the  program  the  non-zero  values  are  stored  in 
an  array  R.  Other  data  read  or  generated  are  the  connectivity 
data  and  material  property  sets.  Connectivity  data  require 
the  element  node  numbers  that  correspond  to  the  nodal  point 
numbers  of  the  complete  assemblage.  Material  property  sets 
identify  different  material  properties  that  can  be 
referenced  by  each  element. 


2.  Element  Stiffness 


element  matrices.  The  nodal  coordinates  for  the  elements, 
material  properties  and  loading  information  which  have  been 
read  and  stored  are  required  here  for  element  matrix 
calculations.  An  element  matrix  is  stored  on  tape  until 
it  is  needed  later  when  assemblage  to  the  structure  matrices 
is  performed. 

For  a  truss  element,  the  element  stiffness  matrix 
is  the  well-known  matrix 


where  A  and  E  are  the  element  section  area  and  E  is  Young's 
modulus  obtained  from  the  material  property  data  input. 

L  is  the  element  length  which  is  calculable  from  the  nodal 
coordinates  and  connectivity  data.  The  element  matrix  is 
expressed  in  terms  of  local  axis.  The  local  element  matrix 
is  transformed  using  direction  cosines  to  obtain  the 
required  global  element  stiffness  matrix  for  use  in  the 
assemblage  of  the  structure  stiffness  matrix. 

3.  Assemblage  of  Structure  Stiffness  Matrix 

The  structure  stiffness  matrix  [K]  is  calculated 
by  direct  addition  of  the  element  stiffness  matrices,  i.e., 

[K]  »  l  [K.] 
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To  perform  the  summation,  each  element  matrix  [K^]  could  be 
written  as  a  matrix  of  the  same  order  as  the  stiffness 
matrix  [K] ,  where  all  entries  in  are  zero,  except  those 
which  correspond  to  an  element  degree  of  freedom.  In 
practice,  only  the  compacted  element  stiffness  matrix, 
which  is  of  order  equal  to  the  number  of  element  degrees 
of  freedom,  together  with  a  connectivity  array  LM  are 
required  to  be  stored.  The  connectiviy  array  relates  to 
each  element  degree  of  freedom,  the  corresponding  assemblage 
degree  of  freedom. 

The  distinct  form  of  the  structural  stiffness  matrix 
lends  itself  to  two  major  reduction  schemes  that  result  in 
saving  computer  storage  space.  One  is  the  storage  of  only 
the  items  within  a  non-zero  band;  and,  by  virtue  of  the 
matrix  being  symmetric,  the  items  in  the  half  band  width  are 
the  only  ones  needed  for  the  analysis.  The  second  is  to 
use  the  "skyline"  technique  whereby  the  columns  above  the 
main  diagonal  are  stored  in  a  one-dimensional  array  of 
length  NWK.  The  column  height  information  is  given  by  array 
MHT.  A  pointer  array  called  MAXA  is  used  to  store  the 
addresses  of  the  diagonal  elements  of  [K]  in  the  A  array. 
Chapter  6  of  Reference  3  and  Chapter  24  of  Reference  7  give 
details  of  these  two  storage  reduction  schemes  which  are 
used  in  SSAP-NPS. 
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B.  SOLUTION  OF  EQUILIBRIUM  EQUATIONS 

Program  SSAP  solves  the  simultaneous  linear  equations 
given  as 


[K]  {U}  =  {R} 

by  using  a  direct  solution  technique  based  on  Gauss 
elimination  (See  Reference  5) .  The  basic  procedure  of  the 
Gauss  elimination  procedure  is  to  reduce  the  coefficient 
matrix  (i.e.,  the  items  in  the  reduced  stiffness  matrix) 
of  the  equation  to  an  upper  triangular  matrix  from  which 
the  unknown  displacements  {U}  can  be  calculated  by  a  back 
substitution.  This  process  is  called  the  Triangular 
Decomposition  of  [K] . 

Subroutine  COLSOL  (see  Section  IV  and  Program  Listing) 

T 

is  used  to  obtain  what  is  called  the  [L] [D] [L  ]  factoriza¬ 
tion  (triangular  decomposition)  of  a  stiffness  matrix  and 
also  to  reduce  and  back-substitute  the  load  vector.  The 
complete  process  gives  the  solution  of  the  element  equili¬ 
brium  equations. 

C.  EVALUATION  OF  ELEMENT  STRESSES 

Once  the  nodal  point  displacements  (U)  have  been 
obtained,  element  stresses  are  calculated  in  the  final 
phase  of  the  analysis.  This  is  accomplished  by  calculating 
the  element  compacted  strain-displacement  transformation 
matrix  and  then  extracting  the  element  nodal  point 
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displacements  from  the  total  displacement  vector  using 
the  LM  array  of  the  element.  By  establishing  the  strain- 
displacement  transformations,  and  relating  stress  to  strain 
by  Hooke's  law,  stresses  can  be  calculated  at  any  desired 
location  in  the  system. 


IV.  PROGRAM  ORGANIZATION 


A.  SUBROUTINES 

Program  SSAP-NPS  is  written  in  an  enhanced  form  of 
BASIC  language  for  use  on  an  H-P  System  45  desktop  com¬ 
puter.  The  program  consists  of  a  main  body  and  14  sub¬ 
routines,  each  of  which  is  designed  to  compute  one  or  at 
most  a  few  basic  steps  in  the  solution  process.  The 
following  is  a  list  of  the  subroutines  and  a  brief 
description  of  their  functions. 

1.  INPUT  -  Reads,  generates,  and  prints 

nodal  point  input  data  to 
calculate  equation  numbers  and 
store  them  in  the  ID  array. 

2.  LOADS  -  Reads  nodal  load  data,  calculates 

load  vector  for  each  load  case 
and  stores  it  on  file  "ILOAD". 

3.  ELCAL  -  Loops  over  all  element  groups 

for  reading,  generating,  plotting 
and  storing  the  element  data. 

4.  ELEMNT  -  Calls  the  appropriate  element 

type  subroutine.  Only  the  truss 
element  is  presently  available. 

5.  TRUSS  -  Sets  up  storage  and  calls  the 

truss  element  subroutine. 

6.  RUSS  -  Truss  element  subroutine  to  read, 

generate  element  information, 
assemble  structure  stiffness  matrix, 
and  performs  stress  calculations. 

7.  COLHT  -  Calculates  column  heights. 

8.  ADDRESS  -  Calculates  addresses  of  diagonal 

elements  in  banded  matrix  whose 
column  heights  are  known. 
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9 .  ASSEM 


10 .  ADDBAN 


11.  COLSOL 


12.  STRESS 


-  Calls  element  subroutines  for 
assemblage  of  the  structure 
stiffness  matrix.  Only  the  truss 
element  is  presently  available. 

-  Assembles  upper  triangular  element 
stiffness  into  compacted  global 
stiffness. 

-  Solves  finite  element  static 
equilibrium  equations  in  core, 
using  compacted  storage  and  column 
reduction  scheme. 

-  Calls  the  element  subroutine  for 
the  calculation  of  stresses. 


13.  WRITE  -  Prints  displacements. 

14 .  PMESH  -  Plots  structure  mesh  and  numbers 

each  node. 


B .  PROGRAM  TAPE 

The  SSAP-NPS  program  is  stored  on  tape  in  two  separate 
files  called  "SSAP"  and  SSAPO".  This  tape,  called  the 
program  tape  hereon,  is  inserted  in  the  H-P  System  45 
standard  tape  drive.  In  using  the  program,  the  "SSAP”  pro¬ 
gram  file  is  loaded  via  the  keyboard  into  the  computer 
memory.  When  the  program  is  run,  it  automatically  gets 
the  "SSAPO"  file  when  it  is  ready  to  use  that  portion  of 
the  program.  Figure  2  shows  how  the  program  segments  are 
organized.  On  the  figure,  the  arrow  indicates  the  part  of 
the  program  to  which  the  "SSAPO"  segment  is  overlayed. 

All  subroutines  below  the  arrow  on  the  "SSAP"  segment  are 
discarded  from  the  computer  memory  since  they  are  no  longer 


needed . 


SSAP  Segment 


MAIN 


ELEMNT 


TRUSS 

RUSS 


"SSAPO"  Segment 


C .  DATA  TAPE 


Before  the  program  can  be  run,  a  data  tape  must  be 
inserted  into  the  H-P  System  45  optional  tape  drive.  The 
data  tape  contains  the  data  files  onto  which  computed  data 
will  be  stored  and  retrieved  by  the  program.  These  files 
are  named  "Al",  MA2",  "A3",  "A4",  "A5",  "ILOAD" ,  " IELMNT" , 
and  "MESH".  The  "MESH"  file  is  required  by  the  plotting 
subroutine  (PMESH) ;  the  other  files  are  discussed  further 
in  the  next  section.  The  data  tape  also  contains  the  files 
which  will  store  input  data  for  the  program.  These  files 
are  named  "COOR",  "LOADS",  and  "TRUSS".  Appendix  A 
discusses  the  input  data  required  for  these  files. 

D.  FLOW  CHART  AND  STORAGE  ALLOCATION 

Figure  3  shows  a  flow  chart  of  the  program.  Figures  4, 
5,  and  6  give  the  data  files  used  and  storage  allocations 
used  during  the  various  program  phases.  Definitions  for 
the  variables  are  given  in  Appendix  A. 


Figure  3.  SSAP-NPS  Flow  Chart 
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Figure  4*  Array  Storage  Allocation  and  Tape  Data 
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Figure  5.  Array  Storage  Allocation  and  Tape  Data 
Location,  Element  Input  Phase 
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Figure  6,  Array  Storage  Allocation  and  Tape  Data 
Location,  Assemblage  and  Solution  Phase 
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V.  PROGRAM  CAPABILITIES 

A.  2-D  TRUSS  STRUCTURES 

The  SSAP-NPS  program  can  solve  for  two-dimensional 
truss  structure  displacements  and  stresses  for  a  given 
loading  condition.  The  program  output  will  first  docu¬ 
ment  the  input  data,  automatically  generating  missing 
input  if  asked  to  do  so  by  the  input  data.  It  will  then 
print  out  a  plot  of  the  structure  in  the  XY  plane,  with 
nodes  properly  labeled  by  node  number.  The  nodal  displace¬ 
ments  in  the  X,  Y,  Z  directions  are  tabulated;  and  finally, 
the  axial  stress  for  each  element  is  tabulated.  A  positive 
stress  denotes  a  truss  element  in  tension;  a  negative 
stress,  an  element  in  compression.  A  2-D  truss  example 
is  included  in  Appendix  B.  The  displacement  and  stress 
calculation  phase  is  repeated  if  more  than  one  loading 
case  is  input. 

B.  3-D  TRUSS  STRUCTURES 

SSAP-NPS  can  solve  three-dimensional  truss  structures 
with  the  same  capabilities  as  for  the  2-D  structures 
discussed  above.  In  addition,  an  isometric  mesh  plot  of 
the  structure  is  outputted.  An  example  in  Appendix  B 
illustrates  a  program  run  for  a  3-D  truss  structure. 
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C .  MESH  PLOTS 


Describing  the  structure  model  or  data  input  seems 
trivial  at  first  but  can  be  the  source  of  much  frustration 
Mesh  points  are  commonly  missed  or  tabulated  incorrectly. 
The  computer  is  unforgiving  and  one  small  error  can  make 
subsequent  work  meaningless.  Some  form  of  data  screening 
is  desirable,  such  as  a  graphical  display  or  plot  of  the 
structure  models  which  can  be  generated  from  a  minimum  of 
geometric  input  parameters.  Mesh  plots  will  verify  that 
the  analyst  and  computer  are,  indeed,  working  on  the  same 
problem.  This  feature  is  only  usable  on  an  H-P  System  45 
computer  that  has  the  graphics  package  option  installed. 

SSAP-NPS  produces  two  types  of  mesh  plot  outputs:  a 
mesh  plot  on  rectangular  XY  axes  for  two-dimensional 
problems  where  all  nodes  are  located  on  the  XY  plane;  and 
an  isometric  mesh  plot  for  three-dimensional  truss 
problems.  The  positive  XY  and  Z  axes  for  the  isometric 
plot  are  120  degrees  apart  with  the  positive  Y  axis  in  the 
vertical  direction.  The  positive  X  and  positive  Z  axes 
are,  respectively,  120  degrees  clockwise  and  120  degrees 
counterclockwise  from  vertical.  In  order  to  avoid  having 
two-dimensional  problems  plotted  on  the  isometric  axes, 
their  nodes  should  be  located  on  the  XY  plane. 
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VI.  REMARKS,  RECOMMENDATIONS  AND  CONCLUSIONS 

A.  RUN  TIMES 

In  using  the  program,  the  size  of  the  data  files  created 
on  the  data  tape  has  a  significant  effect  on  the  amount  of 
time  it  takes  to  complete  a  run.  This  is  because  the  pro¬ 
gram  stores  and  retrieves  data  from  various  locations  on 
the  tape  during  a  run  and  much  time  is  spent  winding  and 
rewinding  the  data  tape. 

As  an  illustrative  case,  take  the  2-D  truss  structure 
example  in  Appendix  B.  With  all  eleven  data  files 
encompassing  approximately  100  defined  redords  on  the 
data  tape,  the  problem  took  about  5  minutes  from  start  to 
finish;  however,  when  the  data  files  encompassed  roughtly 
500  defined  records,  the  problem  took  more  than  9  minutes 
to  finish.  The  program  simply  had  further  to  go  to  access 
the  different  data  files  during  the  execution  of  the  program. 
It  is  therefore  important  to  size  the  data  files  such  that 
they  are  not  much  larger  than  required  to  hold  the  data 
stored  in  them.  This  will  require  some  knowledge  of  mass 
storage  techniques  and  definitions.  One  rough -guide  in 
sizing  the  amount  of  information  a  data  file  can  store  is 
that  a  standard  defined  record  can  hold  256  bytes  of  infor¬ 
mation  and  that  a  real  number  takes  up  8  bytes  while  an 
integer  number  takes  up  4  bytes.  Therefore  creating  a  data 
file  with  5  defined  records  means  that  the  data  file  can 
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n 


hold  approximately  5  *  256/8  real  numbers  or  5  *  256/4 
integer  numbers  or  any  combination  thereof.  For  more 
detailed  information  on  data  storage  and  bytes  per  variable 
type  see  References  4  and  5. 

B .  ERROR  MESSAGES 

From  time  to  time  when  running  the  program,  error 
messages  are  encountered.  When  this  happens,  the  program 
pauses.  In  most  instances,  the  error  can  be  traceable  to 
a  data  input  error,  especially  if  error  17  (subscript  out 
of  range)  is  encountered.  In  cases  when  error  88  (read 
data  error)  or  error  81  (mass  storage  device  failure)  is 
encountered,  it  may  mean  that  the  data  tape  has  been 
stretched  from  continued  use.  In  some  cases  it  is  curable 
by  pressing  the  CONT  key;  otherwise  the  data  tape  may  need 
to  be  re-initialized  or  replaced  and  the  analysis  must 
start  from  the  beginning.  Error  59  (physical  or  logical 
end-of-file)  indicates  that  one  or  more  of  the  data  files 
is  too  small  to  hold  all  of  the  data  to  be  stored  in  it. 
When  this  happens,  the  program  line  number  at  which  the 
error  occurs  is  displayed  by  the  computer  and  the  statement 
will  reveal  the  data  file  for  which  the  error  occurred. 

This  problem  can  be  alleviated  by  enlarging  the  data  file 
size,  i.e.,  purging  the  existing  data  file  and  recreating 
it  with  an  increased  number  of  defined  records.  Other 
software-related  error  numbers  and  their  causes  are 
summarized  in  References  4  and  5. 
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C.  RECOMMENDATIONS 


As  stated  before,  run  times  are  greatly  affected  by 
the  winding  and  rewinding  of  the  data  tape.  One  cure  for 
this  is  to  store  and  retrieve  data  from  disc  mass  storage. 
Using  disc  mass  storage  also  means  that  larger  problems 
may  be  handled  because  of  the  increase  in  storage  capa¬ 
city.  Alterations  to  the  program  should  be  minimal  to  affect 
this  change  in  mass  storage  identifier  for  data  files. 

As  presently  written  the  SSAP-NPS  program  can  only  call 
upon  a  truss  element.  However,  the  program  is  structured 
such  that  other  element  types  can  be  added  in  the  form  of 
additipnal  subroutines.  The  addition  of  other  element  types 
would  also  require  certain  subroutine  modifications  to  allow 
for  the  increase  in  the  nodal  degrees  of  freedom;  and  a  mesh 
plot  subroutine  would  have  to  be  written  for  the  added 
element  types.  The  addition  of  other  element  types  would 
be  a  major  improvement  and  would  make  SSAP-NPS  a  truly 
general  purpose  finite  element  analysis  program.  An  increase 
in  computer  memory  size  would,  in  this  case,  be  imperative 
if  the  program  is  to  be  expected  to  handle  realistically 
large  problems.  This  means  the  addition  of  some  type  of 
disc  mass  storage  medium  and  foregoing  the  use  of  the  data 
tape  altogether. 

The  isometric  mesh  plot  output  is  a  great  visual  aid 
for  an  analyst.  A  plot  of  the  deformed  shape  of  the 
structure  after  the  loads  are  applied  is  another  visual  aid 
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that  would  be  very  helpful  in  quickly  assessing  the  effect 
of  loads  on  a  structure.  In  order  to  effect  this  for  the 
program  as  presently  structured,  it  would  involve  creating 
additional  data  files  and  rewriting  the  plotting  subroutine. 
Another  type  of  plot  that  would  be  helpful  in  quick  analysis 
is  the  stress  plot.  This  would  indicate  at  each  node  loca¬ 
tion  the  magnitude  and  direction  of  the  principal  stress  at 
that  node  so  that  the  analyst  can,  at  a  glance,  locate 
critical  loads  for  a  structure.  The  effort  required  for 
this  improvement  would  be  similar  to  that  in  effecting  a 
deformed  shape  plot. 

D.  CONCLUSIONS 

The  SSAP-NPS  program  is  a  useful  tool  for  finite  element 
analysis  of  truss  structure.  It  has  a  graphics  capability 
for  plotting  an  input  mesh  in  3  orthogonal  planes.  It  also 
has  automatic  mesh  generation  capability. 

Using  a  data  tape  for  storing  and  retrieving  data  is  an 
effective  means  of  reducing  computer  memory  allocation 
requirements;  however,  this  comes  at  a  price  of  increased 
run  times  due  to  tape  winding  and  rewinding  during  a  run. 

It  is  clear  that  in  order  to  make  the  program  more  effi¬ 
cient  that  data  mass  storage  must  be  on  a  disc  storage 
medium. 

The  program  is  structured  so  as  to  make  the  addition  of 
other  element  types  possible  by  appending  new  element 


subroutines.  This  is  a  necessary  step  in  making  the 
program  a  more  general  finite  element  program. 

Additional  graphics  programming  would  also  be  desirable 
in  order  to  fully  take  advantage  of  the  H-P  System  45 
graphics  capabilities.  This  could  take  the  form  of 
effecting  a  deformed  structure  mesh  plot  and  a  stress  plot. 
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APPENDIX  A 


USER'S  MANUAL 


In  this  Appendix  the  data  input  required  for  the 
input  data  files  "COOR",  "LOAD",  and  "TRUSS"  on  the 
data  tape  are  discussed.  The  units  for  any  one  problem 
must  be  consistent.  Error  messages  or  erroneous  results 
in  the  running  of  the  problems  almost  always  indicate  input 
errors  so  careful  preparation  and  review  of  input  data  is 
an  especially  important  task  for  the  user.  At  the  end 
of  the  Appendix  is  an  alphabetical  list  of  variable  names 
and  a  brief  definition  of  each. 

A.  DATA  TAPE  PREPARATION 

Before  the  input  data  can  be  stored,  the  data  tape  must 
have  the  required  data  files  created  on  it,  with  each  file's 
number  of  defined  records  sufficiently  large  to  handle  the 
analysis  at  hand.  Table  I.  is  a  program  listing  which  will 
create  and  size  the  required  data  files.  The  size  indicated 
for  each  data  file  is  strictly  arbitrary  and  were  sized  to 
handle  the  example  problems  in  Appendix  B.  The  user  can 
change  the  record  size  to  suit  a  particular  analysis.  This 
program  must  be  run  with  the  data  tape  in  the  T14  mass 
storage  medium.  If  the  data  tape  already  has  the  data 
files  on  the  data  tape  but  file  sizes  need  to  be  altered. 
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Table  I.  Program  "CREATE"  Listing 


the  program  listed  in  Table  II  will  purge  the  existing 
files  from  the  tape  and  the  file  creation  program  can  then 
be  run  to  recreate  the  data  files  with  the  new  file  sizes. 

B.  DATA  FILE  "COOR" 

This  data  file  contains  the  problem  title,  control 
information  and  the  nodal  point  information.  The  data 
variables  are  recorded  on  the  tape  sequentially  in  the 
following  order: 

HED$,  NUMNP,  NUMEG,  NLCASE ,  MODEX 

N,  ID (1 ,N) ,  ID (2 ,N) ,  ID  ( 3 ,N) ,  X(N),  Y(N),  Z (N) ,  K(N) 

The  second  line  of  data  variables  is  defined  here  as  the 
nodal  information  set  for  node  N. 

1.  Discussion  of  Data  Variables 

HED$  is  the  title  information  and  is  limited  to  50 
alphanumeric  characters,  including  blanks.  The  title  must 
be  enclosed  in  quotation  marks  which  are  not  counted  as 
part  of  the  50  character  limit. 

NUMNP  is  the  total  number  of  nodes  and  controls  the 
number  of  times  the  program  will  look  for  the  nodal  point 
information  indicated  on  the  second  line  of  data. 

NUMEG  is  the  number  of  element  groups.  The  total 
number  of  elements  are  dealt  with  in  element  groups.  An 
element  group  consists  of  a  convenient  collection  of 
elements.  Each  element  group  is  input  as  required  by  data 
file  "TRUSS"  to  be  discussed  later.  There  must  be  at 
least  one  element  group. 
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NLCASE  gives  the  number  of  loading  cases  to  be 
applied  to  the  structure. 

MODEX  is  the  flag  which  indicates  the  solution 
mode  which  gives  the  user  the  option  of  checking  the  data 
without  executing  the  analysis  (MODEX  =  0)  .  If  MODEX  =  1, 
the  execution  mode,  the  program  goes  on  and  solves  the 
problem.  The  program  gives  the  user  the  option  of  over¬ 
riding  the  MODEX  setting  by  pausing  execution  and  asking 
the  user  if  he  would  like  to  continue  the  solution.  This 
is  a  helpful  feature  if  the  user,  after  reviewing  the  input 
data  as  listed  decides  that  the  program  should  go 
ahead  and  solve  the  problem.  The  user  is  thus  saved  the 
trouble  of  changing  the  value  of  MODEX  and  starting  the 
program  from  the  beginning. 

N  is  the  node  (joint)  number  and  must  have  a  value 
1  _<_  N  <_  NUMNP.  Nodal  data  must  be  defined  for  all  NUMNP 
nodes  and  may  be  input  directly  or  automatically  generated 
(see  discussion  of  variable  KN) .  The  last  node  that  is 
input  must  be  NUMNP. 

ID(1,N),  ID(2,N),  and  ID(3,N)  are,  respectively, 
the  X,  Y,  and  Z  translation  boundary  condition  codes.  A 
boundary  condition  code  can  only  have  a  value  of  either  0 
or  1.  A  value  of  0  indicates  unspecified  (free)  displace¬ 
ment;  a  value  of  1  specifies  fixed  displacement,  i.e,,  the 
node  is  not  free  to  move  in  the  specified  direction. 
Concentrated  forces  may  be  applied  only  in  the  directions 
where  ID  has  a  value  of  zero. 


One  system  equilibrium  equation  is  required  for 
each  unspecified  (ID(M,N)=0)  degree  of  freedom  of  the 
structure.  The  total  number  of  equilibrium  equations  is 
defined  as  NEC?  and  is  always  less  than  or  equal  to  3  times 
NUMNP.  Specified  (ID(M,N)=1)  degrees  of  freedom  are  removed 
from  the  final  set  of  equilibrium  equations. 

X (N) ,  Y(N),  Z(N)  are  the  X,  Y,  and  Z  coordinates 
of  node  N. 

The  nodal  information  sets  need  not  be  input  in 
node  order  sequence;  eventually,  however,  all  nodes  from 
N  =  1  to  N  =  NUMNP  must  have  nodal  information  defined. 

Nodal  information  sets  for  a  series  of  nodes 

[NL' N1+1*KN'  N]+2*KN'  •••  #  n21 

may  be  automatically  generated  from  information  given  on 
two  consecutive  sets: 


SET 

1  -  Nx, 

ID (1,N^) , . . . ,  X(N1) , 

.  .  . ,  KN 

SET 

2  -  n2  * 

id(1,n2)  ,  .  .  . ,  x(n2)  , 

.  .  . ,  KN 

KN^  is  the  node  generation  parameter  given  for  the 
first  set.  The  first  generated  node  is  N^+1*KN^;  the 
second  generated  node  is  N^+2*KN2,  etc.  Generation 
continues  until  node  number  N2  -  KN^  is  established.  Note 
that  the  node  difference  N2  -  N^  must  be  evenly  divisible 
by  KN1 . 
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In  the  automatic  generation  of  nodal  information 
sets,  the  boundary  condition  codes  [ID(M,N)]  and  the 
coordinates  [X(N),  Y(N),  Z(N)]  for  the  generated  sets  are 
interpolated  linearly. 

C.  DATA  FILE  "LOADS" 

This  data  file  contains  the  load  data  for  the  structure 
Each  load  case  requires  following  data  entered  sequentially 
The  total  number  of  load  cases  is  defined  by  the  value  of 
NLCASE  previously  recorded  on  data  file  "COOR". 

LL,  NLOAD 

NOD,  IDIRN,  FLOAD 

The  second  line  of  data  is  defined  here  as  a  nodal  load 
data  set. 

LL  is  the  load  case  number  and  must  be  put  in  ascending 
sequence  starting  with  one  and  ending  with  NLCASE. 

NLOAD  is  the  number  of  concentrated  loads  applied  in 
this  load  case  and  defines  the  number  of  nodal  load  data 
sets  that  will  be  input. 

NOD  is  the  node  number  to  which  a  particular  load  is 
applied. 

IDIRN  is  the  degree  of  freedom  number  for  the  load. 

If  IDIRN  =  1,  the  load  acts  in  the  X  direction.  If 
IDIRN  *  2,  the  load  acts  in  the  Y  direction;  IDIRN  =  3, 
in  the  Z  direction. 

FLOAD  is  the  magnitude  of  the  load  and  must  be  acting 
in  the  global  X,  Y,  or  Z  direction  only. 


D.  DATA  FILE  "TRUSS" 

Truss  elements  are  two-node  members  allowed  arbitrary 
orientation  in  the  X,Y,Z  coordinate  system.  The  truss 
transmits  axial  force  only,  and,  in  general  is  a  six  degree 
of  freedom  element  (i.e.,  three  global  translation 
components  at  each  end) .  The  following  sequence  of  data 
is  recorded  on  the  tape  data  file  "TRUSS"  for  each  element. 
The  total  number  of  element  groups  (NUMEG)  was  deifined 
earlier  on  tape  data  file  "COOR". 

NPAR(l),  NPAR(2)  ,  NPAR(3) 

N,  E(N),  AREA (N) 

M,  II,  JJ,  MTYP ,  KG 

The  second  line  is  defined  here  as  a  material/section 
property  data  set.  There  are  NUMMAT  sets  of  material/ 
section  property  data  sets.  The  data  sequence  on  the  third 
line  is  defined  here  as  an  element  data  set. 

NPAR(l)  is  equal  to  one,  the  truss  element  type  number. 
NPAR ( 2 )  is  the  number  of  elements  in  this  element  group. 
NPAR(2)  is  equal  to  the  variable  NUME  in  the  program. 

NPAR (3)  is  equal  to  the  number  of  material/section  property 
data  sets  to  be  input.  NPAR (3)  is  equal  to  the  variable 
NUMMAT  in  the  program,  and  must  have  a  value  equal  to  or 
greater  than  one.  If  a  zero  is  entered,  it  is  defaulted  to 


For  the  material/section  property  data  set,  N  is  the 


number  of  the  data  set  and  must  be  in  ascending  order 

beginning  with  one  and  ending  with  N  =  NUMMAT.  E(N)  is 

Young's  modulus,  and  AREA (N)  is  the  section  area.  The 
Young's  modulus  and  section  area  of  each  truss  element 
defined  by  each  element  data  set  is  identified  by  one  of 
the  data  sets  input  here. 

For  the  element  data  set,  M  is  the  truss  element  number 
and  must  have  a  value  1  <  M  <  NPAR(2),  NUME  elements  must 
be  input  and/or  generated  (see  KN  discussion  below)  in 
ascending  sequence  beginning  with  one. 

II  is  the  node  number  at  one  end  of  the  element.  JJ  is 

the  node  number  at  the  other  end  of  the  element.  II  and  JJ 

must  not  be  equal  and  must  have  values  between  one  and 
NUMNP. 

MTYP  is  the  value  of  material  property  data  set 
applicable  for  the  element. 

KG  is  the  node  generation  increment  used  to  compute 
node  numbers  for  missing  element  data  sets.  If  KN  =  0, 
its  value  is  defaulted  to  KN  *  1.  If  there  are  "J”  missing 
element  data  sets  between  element  M  =  K  and  M  =  L  they  are 
generated  using  MTYP  of  element  number  K  and  by  incrementing 
the  node  numbers  of  successive  elements  with  the  value  KG 
for  element  K.  Table  III  is  a  summary  list  of  variable 
names  used  in  the  program  and  a  brief  definition  for  each. 
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TABLE  III  -  LIST  OF  VARIABLES 


A 

AREA  (N) 

E(N) 

FLOAD 

HED$ 

IDIRN 

ID(P,N) 

II 

JJ 

K 

Ki 

KG 

KN 

LL 

M 

MAXA 

MHT 

MODEX 

MTYPE 


Reduced  stiffness  matrix 

Truss  section  area  for  material/section  property 
data  set  N 

Young's  modulus  for  material/section  property 
data  set  N 

Load  magnitude 

Title  of  problem  ( <_50  characters) 

Degree  of  freedom  number 
If  =  1,  X  direction 
If  -  2,  Y  direction 
If  =  3,  Z  direction 

Translation  boundary  code  for  direction  P, 
node  number  N. 

Node  number  at  one  end  of  truss  element 
Node  number  at  other  end  of  truss  element 
Structure  stiffness  matrix 
Element  stiffness  matrix 

Node  generation  increment  for  element  data  sets 

Node  generation  increment  for  nodal  information 
sets 

Load  case  number  (1  LL  <_  NLCASE) 

Truss  element  number 

Array  to  identify  diagonal  element  addresses 

Column  height  array 

Flag  indicating  solution  mode 
If  =  0,  data  check  only 
If  =  1,  execution 

Material/section  property  set  for  element  M 
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N 


NEQ 

NLCASE 

NLOAD 

NOD 

NPAR(l) 
NPAR ( 2 ) 

NPAR ( 3 ) 

NUME 

NUMEG 

NUMNP 

R 

U 

X(N)  , 
Y(N)  , 


Node  number  (1  <_  N  <_  NUMNP)  or  material/section 
property  data  set  number  (1  <  N  <  NUMMAT) 

Number  of  equations 

Total  number  of  load  cases 

Total  number  of  loads  for  load  case  LL 

Node  number  to i which  load  is  applied 

Truss  element  case  number  (=1) 

Number  of  elements  in  the  element  group; 

NPAR  ( 2 )  =  NUME  >_  1 

Total  number  of  different  material/section 
property  data  set,  NPAR (3)  =  NUMMAT  >_  1 

Number  of  elements  in  the  element  group; 

NUME  >_  1 

Total  number  of  element  groups;  NUMEG  >_  1 
Total  number  of  nodal  points;  NUMNP  1 
Load  vector  array 
Displacements  array 

Coordinates  of  node  N  in  global  system  axes 


APPENDIX  B 


EXAMPLE  PROBLEMS 


In  this  Appendix  the  procedure  in  a  typical  run  for  the 
SSAP-NPS  program  is  discussed.  Then  input  and  output 
for  two  example  problems  are  presented. 

A.  TYPICAL  PROGRAM  RUN 

In  order  to  complete  a  run  for  the  SSAP-NPS  program,  the 
user  must  have  two  tapes:  a  data  tape,  and  a  program  tape. 
The  data  tape  may  or  may  not  have  the  required  data  files 
already  created  on  it.  The  program  tape  must  have  the 
following  programs  on  it:  "SSAP",  "SSAPO",  "CREATE"  "PURGE" 
and  "INPUT".  "SSAP"  is  the  first  segment  of  two  segments 
comprising  the  finite  element  program.  "SSAPO"  is  the 
other  segment  of  the  program.  "CREATE"  is  a  data  files 
creation  program  and  is  listed  on  Table  I  of  Appendix  A. 

The  defined  record  size  for  each  data  file  may  be  changed 
to  suit  a  particular  problem.  "PURGE"  is  a  tape  data  files 
purging  program  and  is  listed  on  Table  II  of  Appendix  A. 

This  program  clears  the  data  tape  of  existing  input  data 
files  if  file  sizes  need  changing.  "INPUT"  is  a  program 
which  writes  the  problem  input  data  on  the  three  input  data 
files  as  discussed  in  Appendix  A.  This  program  varies 
with  each  problem.  Sample  listings  are  given  in  the  coming 
examples . 
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The  following  sequence  of  steps  is  what  is  involved 
in  completing  a  program  run.  Steps  2  and  3  are  omitted  if 
the  data  tape  already  contains  the  data  file  and  the  user 
does  not  want  to  change  any  of  the  file  sizes. 

1.  Insert  the  data  tape  in  T14  mass  storage,  and 
the  program  tape  in  T15  mass  storage.  Turn  on 
computer  power.  Execute  PRINTER  IS  0  statement  if 
Thermal  printer  output  is  desired. 

2.  Load  and  run  "PURGE"  program.  Skip  this  step  if 
the  data  tape  does  not  have  data  files  already 
created  on  it. 

3.  Load  and  run  "CREATE"  program. 

4.  Load  and  run  "INPUT"  program. 

5.  Load  and  run  "SSAP"  program. 

The  program  takes  over  from  this  point.  If  the  value  of 
MODEX  =  0,  the  data  check  only  mode,  the  program  will  pause 
and  give  the  user  the  option  of  either  stopping  the  program 
or  continuing.  If  the  user  continues  the  run,  the  program 
sets  MODEX  =  1  and  completes  the  run.  This  pause  in  the 
program  serves  the  purpose  of  giving  the  user  time  to  review 
input  as  interpreted  by  the  program  and  choosing  either  to 
go  on  with  the  calculations,  or  stop  to  make  corrections 
and/or  changes . 
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B.  3-D  TRUSS  EXAMPLE 

Figure  7  shows  a  three-dimensional  truss  example  with 
loads  in  the  X,  Y  and  Z  directions.  A  tabulation  of  the 
load  magnitudes  acting  on  locations  (nodes)  1,  2,  3  and  4 
are  included  in  the  figure.  Table  IV  lists  the  input 
program  "INPUT"  which  would  be  required  to  accomplish  step  4 
above.  Preceding  comments  are  included  here  only  as  a 
means  to  indicate  which  variables  are  being  defined  by  the 
lines  that  follow.  Table  V  is  the  resulting  output  for 
this  example. 

C.  2-D  TRUSS  EXAMPLE 

The  truss  in  this  example  has  two  load  cases  as  indi¬ 
cated  in  Figure  8.  Table  VI  is  the  "INPUT"  program 
required  for  step  4  as  discussed  earlier.  This  input 
program  illustrates  how  automatic  mesh  generation  can  be 
utilized  to  minimize  input  data.  Preceding  comment  state¬ 
ments  are  included  to  indicate  which  variables  are  being 
defined  by  the  lines  that  follow.  Table  VII  is  the  result¬ 
ing  output  for  this  example  problem. 

Node  locations  for  2-D  trusses  should  be  located  on  the 
XY  plane  in  order  to  avoid  getting  an  isometric  mesh  plot 
in  the  output. 
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Table  V.  Data  Output  for  3-D  Truss  Example 
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Table  VI.  Program  "INPUT"  Listing  for  2-D  Truss  Example 
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2110  MAT  READ  #4;Disp 
2120  FOR  11=1  TO  Numnp 
2130  FOR  1=1  TO  3 
2140  D< I >=0 
2150  NEXT  I 
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